Density matrix renormalization group and wave 
function factorization for nuclei 



in 
O 

o 



in 



(N 
> 
(N 



(N 



T. Papenbrockfl and D. J. Deanf 

f Department of Physics and Astronomy, University of Tennessee, Knoxville TN 
37996-1201, USA 

t Physics Division, Oak Ridge National Laboratory, Oak Ridge, TN 37831, USA 

Abstract. We employ the density matrix renormalization group (DMRG) and 
the wave function factorization method for the numerical solution of large-scale 
nuclear structure problems. The DMRG exhibits an improved convergence for 
problems with realistic interactions due to the implementation of the finite 
algorithm. The wave function factorization of /pg-shell nuclei yields rapidly 
converging approximations that are at the present frontier for large-scale shell- 
model calculations. 
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C3 ' 1. Introduction 



The nuclear shell model with full configuration mixing is the adequate tool for a 
quantitative description of light- and medium-mass nuclei. In recent years, the no- 
core shell model has successfully been applied to investigate p-shell nuclei |2I, 
f~| , and full-space diagonalizations involving model spaces with up to one-billion Slater 

determinants are now possible for sc?-shell and p/-shell nuclei within the traditional 
shell model Note that it took more than a decade to go from sd-shell nuclei to 

' p/-shell nuclei, and progress was due to more sophisticated algorithms and an increase 

^ , in computer power. The description of heavier nuclei or drip-line nuclei, however, is 

based on increasingly larger model spaces. For such nuclei, exact diagonalizations will 
be unavailable for the next years, and one has to introduce approximations. 

Several authors have suggested truncations of the model space. In many of these 
methods, the selection of the relevant basis states is based on physical insights and 
arguments and is done "by hand" (HIIZIIHI- In other methods, the Hamiltonian itself 
selects the most important basis states. One example is the Monte Carlo shell model 
in], where the huge Hilbert space is sampled stochastically, and only relevant basis 
states are kept. Another example is provided by the density matrix renormalization 
group (DMRG) |10[ [TT] and the wave function factorization j^, where the most 
important basis states are determined from a variational principle. It is the purpose 
of this article to describe recent progress with the latter two methods. We report on 
DMRG results for sd-shell and ^i/-shell nuclei, and use the wave function factorization 
to compute accurate approximations for low- lying states in O/5/2 shell nuclei 

that are at the frontier of full space diagonalizations. 
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2. Density matrix renormalization group 

During the last decade, the DMRG has become the method of choice to compute 
ground state properties for spin-chains. We refer the reader to the review ^Hl and 
briefly describe the main aspects of this method. In the one-dimensional spin systems, 
the lattice is divided into two parts (e.g. "left" and "right"), and the ground state 
energy is computed for a small number of lattice sites. The ground state can be 
expressed in terms of basis states |^) and \r) of the "left" and "right" part of the 
chain, respectively, as 



In the DMRG, one seeks an approximation of j'l') in terms oi m < M basis states. 
These states are the eigenstates of the density matrices pw — Ylr '^ir'^i'r and 
Prr' — X); ^ir^'ir' that corrcspond to the m largest eigenvalues. This is a truncation 
for which the difference between exact ground state and approximated ground state 
has minimal norm One then retains these m states for the left and right part 

of the lattice, and adds another lattice site to either part. This procedure defines 
the "infinite algorithm" and is repeated until the desired system size is reached. At 
this point, an iterative procedure known as the "finite algorithm" is used to enlarge 
one part of the lattice at the expense of the other while keeping the total number of 
lattice sites constant. This defines a "sweep" through the system. Typically one finds 
that the truncation error decreases exponentially with the number m of kept states. 
Note that one does not store the wave function in the DMRG. Instead, one keeps 
m-dimensional matrix representations of those operators that define the Hamiltonian 
and other observables of interest. 

Only recently has the DMRG also been applied to finite Fermi systems, and we 
refer the reader to the recent review by Dukelsky and Pittel 14 . There are at least 
three novelties for the DMRG in finite Fermi systems. First, an order of lattice sites, 
i.e. single-particle orbitals has to be chosen. Second, a partition of the system in 
two parts has to be chosen. Third, the two-body interaction induces interactions 
between very distant single-particle orbitals and thereby differs from the spin chains 
with neighboring interactions. For nuclei, Dukelsky and coworkers suggest an ordering 
of the single-particle orbitals based on their distance from the Fermi surface and 
based the partition on particle versus hole orbitals 15 . Using the infinite DMRG 
algorithm, they obtained rapidly converging results for schematic pairing and pairing- 
plus-quadrupole models but encountered rather slow convergence for sc?-shell nuclei 
with more realistic (and complex) interaction |15| . 

Our implementation of the DMRG is as follows. We choose a partition of neutron 
orbitals versus proton orbitals. This is motivated by our experience with the wave 
function factorization. The spherical single particle orbitals have quantum numbers 
a = (n,l, jz,Tz). We work in the m scheme and conserve total angular momentum 
Jz at each step of the DMRG. As the equivalent of a single lattice site, we choose 
conjugate pairs (a, a) of single-particle orbitals. This is designed to improve angular 
momentum properties and to properly treat pairing correlations. The conjugate pairs 
of single-particle orbitals are ordered such that the most active orbitals, i.e. the 
orbitals close to the Fermi surface are in the center of the chain of orbitals. This 
is motivated by the extensive studies reported by Legeza on this subject To 
overcome difficulties with the long-ranged two-body interaction, we use the finite 
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algorithm. We do not use the infinite algorithm. Instead, we start the finite algorithm 
based on the spherical Hartree-Fock configuration and its lp-\h excitations. We find 
that this approach is superior to the infinite algorithm and the approach suggested by 
Xiang for treating systems with conserved quantities ^Sl- Note that we also considered 
DMRG calculations in a Hartree-Fock basis. Compared to the DMRG in the spherical 
basis, this approach lowered the energy at the start of the finite algorithm. However, 
the breaking of axial symmetry led to larger-dimensional matrix problems and did not 
improve the convergence of the method. In the practical calculations, we start with 
the finite algorithm with a rather small number m of kept states and increase this 
number every two sweeps. 

We performed DMRG calculations for the sd-shell nucleus ^*Si using the USD 
interaction ^H]- The order of the (pairs) of single-particle orbitals is ^38, dsl, ^56, 
dsS, sil, ^5!, ^5!, sil, ^58, dsS, dsl, ^38, where we used the notation l2j\2jz\, and it is 
understood that the first half of orbitals corresponds to proton states, while the (mirror 
symmetric) second half of orbitals corresponds to the neutron states. The initial set 
of states consists of the (^5/2)^^ configuration and its lp-\h excitations. The left part 
of Fig. n shows the DMRG results for the ground state energy as a function of the 
number of states we keep. The energy converges exponentially rapidly as more states 
are kept, and the true ground state energy can accurately be determined from an 
exponential fit to the data. The inset shows the dimension of the DMRG eigenvalue 
problem as a function of the number of kept states. Note that the full eigenvalue 
problem has dimension d = 98710. 




Figure 1. Left: DMRG ground state energies for ^*Si as a function of tlie 
number m of kept states (data points connected by dots), exact result (daslied 
line), and exponential fit (full line). The inset shows the dimension d of the 
DMRG eigenvalue problem vs. m. Right: Same as left except for ^^Ni. 



This promising result motivates us to perform DMRG calculations for the /p-shell 
nucleus ^^Ni. We use the KB8 interaction ,20. .21} . The exact diagonalization of this 
shell-model problem has to deal with about one billion Slater determinants and has 
only been accomplished recently The order of single-particle orbitals is /sS, 758, 
/5I, Pil, /y^, /tS, P33, /yS, /7I, /7I, /tS, P^i, /yS, AT, Pil, /5I, /sS, /sS. 
The initial set of states consists of the {fj/2)^^ configuration and its Ip-lh excitations. 
The right part of Fig. ^ shows the results. The DMRG results seem to be converged 
as an increase of m does not significantly lower the energy. However, the true ground 
state energy is almost 400 keV lower in energy. We recall that the structure of ^^Ni 
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is to some extent similar to ^^Si. This similarity has motivated our choice of states 
to be kept at the beginning of the finite algorithm and our order of the single-particle 
orbitals. It is thus unexpected that the DMRG converges for ^*Si but fails to fully 
converge for the computationally larger problem ^^Ni. 

Note however that our DMRG computations exhibit an improved convergence 
when compared to the pioneering study [Ifij . The most important new ingredients 
are the use of the finite algorithm, the abandoning of the infinite algorithm, and the 
ordering of the single-particle orbitals. We expect that further improvement of the 
convergence relies on an optimal order of single-particle orbitals and on optimal initial 
conditions for the finite algorithm. The method presented in the following section 
avoids these issues and makes use only of the unproblematic truncation step of the 
DMRG. 



3. Wave function factorization 



Modern shell model codes build their basis from products of Slater determinants [tt) 
for the proton space and Slater determinants \v) for the neutron space. Note that the 
dimensions of the proton and neutron space are modest compared to the dimension 
of the product space. A shell-model ground state might thus be expanded as 

|vl')-5]*,.|7r)|^), (2) 

where the sum runs over all available proton and neutron determinants. This 
expansion is similar to Eq. (QJ used in the DMRG. The amplitude matrix ^f^i, is 
rectangular in general and can be factored as — U SV'^ by means of the singular 
value decomposition (SVD). Here, U and V are unitary matrices that operate on the 
proton space and the neutron space, respectively, and S is zero except for nonnegative 
entries Sjj > (the "singular values") along its diagonal. The SVD thus yields 
correlated proton and neutron states 

\P,) = Y.^.j\7r), (3) 

TT 
TT 

and the ground state becomes |^) = J^j ^j\Pj)\^j)- Note that the DMRG is closely 
related to the SVD: ^E'^E't = USSW' and ^ftvl/ = VS^'SV'' are density matrices which 
are diagonalized by the matrices U and V, respectively, and their eigenvalues are 
the squares of the singular values JI]. Note also that the SVD yields exponentially 
rapidly decreasing singular values for ground states obtained from realistic shell-model 
interactions [52|. This observation is interesting for two reasons. First, it shows that 
the DMRG truncation step also works for nuclei. Second, it provides the basis for 
a powerful approximation and basis state selection scheme. In the wave function 
factorization, we approximate the ground state as a sum over products of correlated 
proton and neutron states 
n 
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Note that the states \pj) and \nj) are not normahzed. Variation of the energy yields 
the foUowing set of eigenvalue equations that determine the states \pj) and \nj) 
n 



These equations can be solved iteratively. One starts with a random set of orthogonal 
neutron states \nj) and solves the first of the eigenvalue problems in Eq. lO for 
those proton states that correspond to the lowest energy E. These are then inserted 
in the second eigenvalue problem in Eq. ((HJ, which is solved for the neutron states 
corresponding to the lowest energy. This procedure is repeated until the energy E 
converges. For even-even nuclei, convergence is typically reached in 5-10 iterations, 
while other nuclei require a somewhat larger number of iterations. The eigenvalue 
problems have dimension QDp and fi-Dn for the proton states and neutron states, 
respectively, where Dp (Dn) denotes the dimension of the proton (neutron) space. 
Typically, one needs only a few hundred states for more than 99% overlap with 
the exact ground state. Thus, <C Dp, Dm and the factorization requires one to 
solve eigenvalue problems of small dimensions relative to the exact diagonalization 
which has a dimension that scales as DpDn- Note that the proton and neutron states 
that enter the eigenvalue problem © can be kept orthogonal, and this reduces the 
generalized eigenvalue problem to a standard eigenvalue problem. Note also that axial 
or rotational symmetry can be conserved within the ansatz For details, we refer 
the reader to Ref. 22 . 

In earlier work, we applied the factorization to sc?-shell nuclei and /p-shell 
nuclei and demonstrated its accuracy and capability by comparing with results 
from exact diagonalization. In this work we use the factorization and perform 
structure calculations for a few A — 76, 78 nuclei. The model space consists of 
the 0/5/2 ipOgg/2 orbitals, and the interaction is a monopole-corrected G-matrix [21^ . 
taking 28^128 as a closed core. These problems are at the present frontier of nuclear 
structure calculations. We mention, for instance, the recently accomplished exact 
diagonalziation for double beta-decay studies in ''^Ge and ""^Se |24[ I25j . and the 
determination of the T = 1 two-body matrix elements by fit |26| . 

In this article, we focus on the nuclei ''^Ge, ^^Se, and ^^Kr, and compute their 
lowest-lying states using the m-scheme factorization. For Se and Kr, we also use 
parity to reduce the dimensionality of the shell-model problem. Figure [3 shows the 
energy difference E{U) — Esm of the ground states as a function of the number 
of kept states. The shell- model energies Esm are determined by an exponential fit 
to the numerical data points and are listed in Table [H Note that the wave function 
factorization uses a few hundred states out of up to several hundred thousand proton 
and neutron states. Deformed nuclei typically require more states as deformation 
is driven by proton-neutron interactions j27|. and the corresponding wave functions 
exhibit stronger entanglement between proton states and neutron states. Nevertheless, 
the wave function factorization is very efficient. For ^*Kr, for instance, the model 
space contains about 320,000 Slater determinants for the protons and neutrons, and 
these combine to about a 1.4 x 10^ dimensional configuration space of definite parity. 
The largest eigenvalue problem we solve within the factorization uses only about 400 
correlated states for protons and neutrons, has a much smaller dimension of 3.5 x 10^, 
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Figure 2. Data points: Ground state energies as a function of the number of 
retained states. Dashed lines: exponential fit to data. 



and the Hamiltonian matrix contains about 2.6 x 10^ nonzero matrix elements. These 
numbers are small compared to what exact diagonalizations for comparably large shell 
problems would require [23 • Note also that the angular momentum expectation value 
( J^) is practically an integer value for the largest numbers of retained states. 

We compute excited states as a by-product of the ground state factorization. This 
approach yields the optimal proton (neutron) configurations of the excited states in 
the presence of the neutron (proton) configuration of the ground state. The excitation 
energies of the first excited J'^ = 2+ states are also listed in Table They are 
obtained from the largest calculations we performed and are expected to be upper 
bounds (as the ground state converges more rapidly than the excited states), and 
the theoretical uncertainty is about 150 keV. The comparison with the experimental 
results is reasonable, and for higher precision of the theoretical results one would need 
to target excited states directly [32]. 



Table 1. Shell-model energy -EgM of the ground state and the excitation energies 
£'2+ of lowest J'^ = 2+ state compared to experimental results E'^ . Unit of 
energy is MeV. 



Nucleus 


EsM 






76 Ge 


-57.23 


0.70 


0.563 


76 Se 


-75.74 


0.46 


0.559 


78Kr 


-96.06 


0.44 


0.455 
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4. Summary 

In this article we described recent progress with the DMRG and the wave function 
factorization. Within the DMRG, we obtained fairly well-converged results for large 
scale nuclear structure problems using realistic interactions. Keeping 100-200 states 
yielded energy deviations of a few hundred keV. The improved convergence rests 
mainly on the implementation of the finite DMRG algorithm, an improved ordering 
of single-particle orbitals, and the use of well-selected states at the beginning of the 
finite algorithm. Further progress should be possible by optimizing these aspects. 

We used the wave function factorization for nuclear structure studies of a few 
0/5/2 Ip 0(79/2 shell nuclei. Typically, one needs only a few hundred correlated 
proton and neutron states for an accurate approximation of low-lying states, and 
the error involved is of the order of 100 keV. Our results show that accurate structure 
calculations for fpg-sheW nuclei and tests of the effective interactions are possible at 
a fraction of the effort associated with exact diagonalizations. 
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